Front Propagation in Reaction- Superdiffusion Dynamics - Taming Levy Flights with Fluctuations 
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We investigate front propagation in a reacting particle system in which particles perform scale-free random 
walks known as Levy flights. The system is described by a fractional generalization of a reaction-diffusion 
equation. We focus on the effects of fluctuations caused by a finite number of particles. We show that, inspite 
of superdiffusive particle dispersion and contrary to mean field theoretical predictions, wave fronts propagate at 
constant velocities, even for very large particle numbers. We show that the asymptotic velocity scales with the 
particle number and obtain the scaling exponent. 
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One of the fundamental processes involved in non- 
equilibrium pattern formation is the spatial propagation of in- 
terfaces or fronts. Front propagation usually emerges when 
a local reaction dynamics interplays with diffusion in space 
of the reacting agents and has been observed in a wide range 
of physical, chemical and biological systems fl*]. Prominent 
examples are patterns of bacterial colonies |2], spreading phe- 
nomena in population genetics (sl 0] and interface dynamics 
in superconducters 1 5 ] . One of the most prominent models 
which displays propagating fronts is the Fisher-Kolmogorov- 
Petrovsky-Piscounov (FKPP) equation for the spatial concen- 
tration u{x,t) of a reacting agent, 

dfU = '-fu{l — u) + Au, (1) 

where diffusive motion of the reacting agents is assumed. 
However, this assumption cannot be justified for a number of 
systems. In fact, superdiffusive dispersion in space has been 
observed in a wide range of physical and biological systems, 
e.g. intermittent chaotic systems |6], bacterial motion jTU, and 
foraging patterns of albatrosses |8]. 

Superdiffusive stochastic motion is usually characterized 
by a lack of scale in the microscopic single step distribution. 
One of the most successful theoretical concepts devised for 
the understanding of superdiffusion is a class of random walks 
known as Levy flights ^ . A Levy flight consists of random 
single steps Ax which are drawn from an inverse power-law 
pdf (probability density function) p( Ax) ^ |Aa;|~'^+'^' char- 
acterized by a Levy exponent < /i < 2. Due to the heavy 
tail, the variance in step size is divergent, the process lacks 
a spatial scale and the position X{t) of a Levy flight scales 
heuristically with time t as X{t) t^/^. The associated 
diffusion equations contain fractional generalizations of or- 
dinary derivatives I lol [Till . These fractional Fokker-Planck 
equations can exhibit behaviors strikingly different from or- 
dinary ones 1 12] and have found wide application in physics, 
e.g. protein motion on folded hetero-polymers 1 13] and the 
dynamics of modern epidemics Il4ll5ll . 

In two recent studies wave front dynamics was shown to be 
drastically different from ordinary reaction-diffusion dynam- 
ics when the reacting agents move superdiffusively 
The authors considered fractional generalizations of Eq. Q 



and showed that the spatio-temporal shape u{x, t) of the lead- 
ing egde of a propagating front has a power-law tail along 
the spatial coordinate and accelerates exponentially in time 
as opposed to the constant velocity and exponential decay in 
space exhibited by ordinarily diffusive system. The predic- 
tions made by mean field theory thus indicate that scale-free, 
superdiffusive dispersion of the reacting agents excludes con- 
stant velocity wave fronts and induces an entirely different 
spatio-temporal behavior. 

However, as has been shown in a number of recent studies, 
the effect of fluctuations can be rather profound in these sys- 
tems (Tsi fiol Eoll . For instance, fluctuations can destabilize 
homogeneous wave fronts in two-dimensional systems |21]. 
A finite albeit large number N of particles or reacting agents 
leads to a multiplicative noise term in the reaction-diffusion 
equations and although the variance of the noise is of order 
1/iV, corrections to macroscopic quantities such as the front 
speed V scale as (logA^)^^, a significant correction to the 
mean field approximation even for very large values of A^. 
Brunei and Derrida i23.E3ll extended mean field dynamics by 
an effective cutoff parameter e for the concentration of parti- 
cles below which no reaction and hence no exponential growth 
of the leading edge of a front is possible. Despite the fact that 
a rigorous equivalence with multiplicative noise is still lack- 
ing, the effective cutoff approach is very intuitive and in re- 
markable agreement with simulations of the full probabilistic 
dynamics. 

Here, we focus on the effect of fluctuations on reaction- 
superdiffusion kinetics. We show that for arbitrarily small 
fluctuations (i.e. arbitrarily large particle numbers N), wave 
fronts propagate asymptotically at constant velocities. Fur- 
thermore, we show that as soon as fluctuations enter the de- 
scription the algebraic tail along the spatial coordinate of the 
leading edge disappears and is replaced by an exponential de- 
cay. Thus, despite the fact that reacting agents move superdif- 
fusively in space, the wave front patterns are qualitatively the 
same as in the ordinary diffusion case. We show that a front 
speed V is selected after a transient time and that v scales with 
particle number as w N'^/t^ for Levy exponents fi < 2. 
The results reported here are rather counterintuitive, deviate 
strongly from the predictions cast by mean field theory, and 
indicate that fluctuations affect reactions-superdiffusion sys- 
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Figure 1 : Propagation of wave fronts in superdiffusive systems witli 
effective cutoff, (a) Location of the wave front x*{t) as a function 
of time t for different particle numbers TV (solid lines) for jj. — 1.5 
and 7 — 0.1 in Eq. <11> . Following initial transients the velocity 
V — dx*{t)/dt is constant (dashed lines) and increases with A'^. The 
inset depicts a magnification of the initial phase. The wave front 
accelerates exponentially but deviates from the mean-field dynamics 
(e = 0, dashed line) after the transient period, (b) Shape of the wave 
front at exponentially increasing time steps t = 1.5™ with m — 
5, 6, 20 for A'^ — 10^. After a transient phase the shape remains 
unaltered, decays sharply for 1/A^ < u{x,t) < 1 and follows a 
power law u{x, t) ~ a;-<i+'') for large x. The dashed line indicates 
the effective cutoff e — 10~*. The spatial extent in the numerical 
integration was L = 2^^ f» 8.38 x 10*^. 



tems severely and may by no means be neglected. 

We begin with a simple two particle type { A, B) reaction 
scheme. 



ki 



A^ + B^ 



fc2 



Ax , Bx > Ay , By . 



(2) 
(3) 



Particles of type A and B either react to produce two parti- 
cles of type A or two particles of type B at rate fci and k2, 
respectively. Furthermore, particles of both types may jump 
from position x to position y with a probability density rate 
/ ( I a; — y\) which we assume to be a decreasing function of 
distance \x — y\. The dynamic stochastic variables are the par- 
ticle numbers nA{x, t) and nB{x, t) of particles in an volume 



of size around x of type A and B, respectively. The volume 
il is assumed to be large enough to contain a large number 
of particles but small compared to the spatial extend of the 
system. The total number of particles in Vl around x is given 
by N{x,t) = nA{x,t) + nB{x,t). Without spatial disper- 
sion of particles the local dynamics is governed by the master 
equation 

dtp{n,t) — w{n\m)p{ni, t) — w{m\n)p{n^t) (4) 

m 

for the probability p(n, t) of finding a number n = ua 
of particles of type A at a location x with initial condition 

p(n, 0) = 5{n — no) and the rate 

w{n\m) = ^m{N-m)5n,m+i + ^m{N~m)Sn,m~i- (5) 

The dynamics of the expectation value of {n{t)) is governed 
by 

dt {n{t)) = ^ {n{t) {N - n{t))) , (6) 

where 7 ~ ki — k2- Note that {n{t)) is continuous in [0, N]. In 
a spatially extended system the number of particles is a func- 
tion of position, i.e. n = n{x, t). Apart from normalization, 
(ri(a;, <)} may be interpreted as the probability of finding an A 
particle in the volume $7. Dispersion contributes to the change 
of particles dt {n{x)) according to 

dt{n{x,t))^ f dyfi\x-y\)[{niy,t))-{nix,t))], (7) 



incorporating the probability density rate (Eq. (|3}) of jump- 
ing from y to x. Eq. (0 defines the operator C acting on the 
field {n{x,t)). Denoting the spatial density of particles by 
u{x, t) = n{x, t)/n and combining Eq. with |7) yields 



dt {u) ^-f{u{l-u))+C(u) 



(8) 



where we have set without loss of generally the maximum 
density N/ ft to unity. 

When particles perform ordinary random walks, i.e. if the 
jump rate f{\x — y\) is equipped with a length scale cr, the 
operator C can be approximated by the ordinary Laplacian A 
on scales larger than a. However, this description is no longer 
valid if f{x) ^ |a;|"(^+^' with < ^ < 2. In this case 
individual jumps lack a scale, particle perform Levy flights, 
and C is proportional to a non-local singular integral operator, 



A>^/^uix,t) ^ / dy 



[u{y,t) - u{x,t)] 
I.T-7;|i+A' 



(9) 



where (7^ is a constant |24]. In Fourier space the operator 
^fj./2 equivalent to a multiplication by |fc|^, generalizing the 
well known factor corresponding to the ordinary Laplacian 
which is why A'^/^ is frequently referred to as a fractional 
Laplacian. 
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Figure 2: Scaling beiiavior of tlie fronts profiles at exponentially in- 
creasing times t as in Fig.^in the reference frame of the front. Ar- 
rows indicate the temporal direction. The curves approach a steady 
state characterized by an exponential decay in space, exp — Aa::. The 
position of the from x*{t) is given by u[x*{t),t) = 0.05 . Other 
parameters are as in Fig.Q 

Since (v?^ > (it)^ one cannot computed the dynamics 
of {n{x,t)) from Eq. (|8} alone. The crudest approximation 
is the mean field approach in which fluctuations are ignored 
i{u^) ~ {M)^and thus u = (u)) and Eq. (|8} reads 

5tu==7u(l-it) + A'^/^u, (10) 

which is the Fisher-equation Q for fi — 2 and its superdiffu- 
sive mean field generalization for /i < 2. 

The next better step is to account for fluctuations (A^ 3> 1) 
which leads to an additional multiplicative white noise term 
of variance u(l — u)/N in Eq. MQ\ . Alternatively one may 
incorporate the finite particle number as an effective cutoff 
e ~ 1 /N in the logistic growth term retaining a deterministic 
equation, i.e. 

5tu = 7u(l-M)e(u-e) + A''/2u, (11) 

in which 8 is the Heaviside-function. Qualitatively, this ac- 
counts for the fact that no growth on average can occur if the 
particle concentrations u is less than one particle per unit vol- 
ume. This approach has been applied successfully in the or- 
dinary diffusion scenario |22]. In Eq. (Ill> . mean field theory 
implies e = 0. 

In the following we investigate front dynamics of the sys- 
tem defined by Eq. ( II H and compare its properties to the mean 
field dynamics (e = 0) , the limiting case of ordinary diffusion 
{pL = 2) and the dynamics of the full probabilistic dynamics 
(Eqs. and Q). 

Fig.^depicts front dynamics obtained from numerical in- 
tegration of ( II 1> for a concentration t) initially peaked at 
the origin. For various particle numbers N the velocity and 
shape of the front were computed. Fig.fQa) shows that even 
for very large particle numbers wave fronts move asymptot- 
ically at constant speeds, a remarkable difference from the 
mean field limit (dashed line in the inset) which predicts ex- 
ponentially accelerating front. Fig. [Qb) depicts snapshots of 



the wavefronts on a double-logarithmic scale. After a tran- 
sient phase, the shape of the front approaches a steady state 
with a sharply decreasing boundary at intermediate values for 
the concentration and an algebraic tail for large x with con- 
centrations u{x,t) below the cutoff e. Fig. |2] displays the 
steady state front shape on a semilogarithmic scale. As op- 
posed to the algebraic tail predicted by mean field theory, 
the boundary is an exponential function of relative position 
z — X — vt, i.e. u{z) cxp— Az, characterized by a spa- 
tial scale A. In summary, the characteristic spatio-temporal 
wave front solution of Eq. (II 1> for large times is given by 
u{x, t) ^ exp [— A(x — vt)] , for 1/iV < w < 1 followed by a 
power-law tail u{x,t) ~ (x — for u < 1/N. The 

decay parameter A and the velocity v depend on the particle 
number N. Qualitatively, the this dependence can be deter- 
mined in the moving reference frame under the assumption 
that u{x,t) ~ u{x — vt) — u{z). Figs.fQb) and |2] suggest 
that u{z) = ui{z) = exp — Az for z < z* and e < u ^ 1 
and u{z) = 1*2(2) = A/z^^'^ for z > z* {u < e) where z* 
marks the crossover between exponential and algebraic decay. 
Inserted into Eq. (II 1> yields 

vXui « ui + A^/^M z < z* (12) 
-vu'^ = A^'/^M z > z*. (13) 

Near the crossover the values of the dispersion in both equa- 
tions are approximately the same, i.e. (A''/2w)(z* - Az) « 
(A''/^u)(2* + Az) and both equations can be combined to 
vXui K, ui ~ VU2- With ^2(2*) = —(1 + fi)u2{z*) / z* and 
U2{z*) « wi(z*) we obtain vA « 1 — v{l + /i)/z*. Since 
z* ^ 1 the second term can be neglected and we find 

?;-A"^ (14) 

The dependence of A on can be obtained by the short time 
dynamics. Consider an initial condition u{x, 0) = exp(— Ax) 
on the half-line a; > 0. After a short time Ai ^ 1 the 
wave front is approximately given by u(a;, At) w exp — Aa; + 
At/Xx^~^^. The crossover x can be defined as the point at 
which both terms are of the same order of magnitude, implic- 
itly provided by Acxp— Aa; = At/a;^+'^. Since A <C 1 one 
can approximately solve for the crossover, x « ln(A^ At) / A. 
In order for the exponential to remain invariant under the dy- 
namics, A must be chosen such that the crossover coincides 
with the effective cutoff, i.e. u{x, At) = e which implies 
the scaling relation A (x e^^'^. With Eq. il4\ one obtains the 
scaling law 

wociV^/'^ (15) 

for the velocity v. 

Fig- El shows the front velocity v as a function of particle 
number N obtained by numerical integration of the dynamics 
(Eq. (Ell)- The numerics agree well with the scaling law dlSI l 
over several orders of magnitude and several choices of the 
Levy exponent /i. 
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Figure 3: Asymptotic front velocity w as a function of the particle 
number for different Levy exponents = 1.2, 1.5 and 1.8 (cir- 
cles, diamonds and triangles, respectively). The dashed line indicates 
the scaling v{N) ~ TV^^^. 



the front velocity v is constant after a transient phase, the ve- 
locity obeys the scaling law M5\ and the wave front decays 
exponentially in space. 

We are convinced that our results are of major importance 
for the understanding of front propagation in pattern forming 
systems in which the reactive agents defy the rules of ordi- 
nary diffusion. We have shown that constant velocity fronts 
are typical for pulled front dynamics, contrary to what is ex- 
pected from mean field approximations and we believe that 
our results will contribute to the understanding of more com- 
plex pattern forming systems such as the geographic spread of 
human epidemics 1151. 
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Figure 4: Asymptotic front velocity in the full stochastic model. The 
solid line depicts the total mass I{t) type A particles as a function of 
time. The system size is L = 10*, the total number of particles in 
the system 2.5 x 10''. Other parameters are 7=1 and ^ = 1.8. A 
constant front velocity (constant AI(t) /Af) is attained asymptotically 
(dashed line). The upper inset depicts the front velocity v{N) as 
a function of particles per site for A'^ = 10, 25, 50, 100, 200 and 
1000 and the scaling law <15t (dashed line). The lower inset depicts 
an average n{x) of A particles at position x at an intermediate time t. 
The average was computed from a sample of 100 realizations of the 
process. The dashed line indicates an exponential decay, the dotted 
line the average concentration of one particle per site. 

The front propagation characteristics of the macroscopic 
description in terms of the reaction-superdiffusion equation 
with effective cutoff (Eq. ( II \V ) coincides with our simulations 
of the microscopic dynamics defined by Eqs. (|2j and (|3j- The 
simulation results are summarized in Fig. |3 As in Fig. ^a) 
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